#------------------------------------------------------------------------------

rm(list=ls())
library(LalRUtils)
libreq(tidyverse, data.table, 
       zoo, tictoc, fst, fixest, 
       PanelMatch, patchwork,
       rio, magrittr, janitor, did, 
       panelView, ggiplot, 
       tictoc, binsreg, interflex)
set.seed(42)
theme_set(lal_plot_theme())
#install.packages("ggiplot", repos = "https://grantmcdermott.r-universe.dev")

#------------------------------------------------------------------------------



#------------------------------------------------------------------------------
# Define Paths
#------------------------------------------------------------------------------

# R studio
setwd( dirname(rstudioapi::getActiveDocumentContext()$path) )
# R default : unccoment if you use default R
# setwd(getSrcDirectory(function(){})[1])
#------------------------------------------------------------------------------


#------------------------------------------------------------------------------
# Load Data
#------------------------------------------------------------------------------

vcf <- fread("vcf_data_complete.csv", sep = ",")


#------------------------------------------------------------------------------


#------------------------------------------------------------------------------

# %%
if (exists("vcf_data")){
  vcf = vcf_data[year>= 1995]
}
vcf[, t := year - 1995]
vcf[, time := year - first_pesa_exposure]


# %% above median sample - cutoff in 1990
evsamp = vcf[cover_1990 > quantile(vcf$cover_1990, 0.5)][time %between% c(-6, 10)]
es1 = feols(forest_index  ~ i(time, sch,  ref=-1) | cellid + styear,
            cluster = ~ blk, evsamp)
es2 = feols(green_index ~ i(time, sch,  ref=-1) | cellid + styear,
            cluster = ~ blk, evsamp)
es3 = feols(built_index ~ i(time, sch,  ref=-1) | cellid + styear,
            cluster = ~ blk, evsamp)
# %%
ggi = function(x) ggiplot(x, theme = lal_plot_theme()) + labs(y = "", x = "Time")
(ff = ((ggi(es1) + ggtitle("VCF Forest Index") +
         theme( legend.position = 'none', 
                axis.title = element_text( size = 22 ), 
                axis.text = element_text(size = 20) ) ) /
        ( ggi(es2) + ggtitle("VCF Non-Forest Green Index")  +
         theme( legend.position = 'none', 
                axis.title = element_text( size = 22 ), 
                axis.text = element_text(size = 20)) ) /
         (ggi(es3) + ggtitle("VCF Bare Ground Index") +
         theme( legend.position = 'none', 
                axis.title = element_text( size = 22 ), 
                axis.text = element_text(size = 20)))
         )
)
ggsave("appendix_figurea4.pdf", ff, width = 7, height = 10,
       device = cairo_pdf)


#------------------------------------------------------------------------------

